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We investigate the spectroscopy and decays of the charmonium and upsilon systems in a potential 
model consisting of a relativistic kinetic energy term, a linear confining term including its scalar 
and vector relativistic corrections and the complete perturbative one-loop quantum chromodynamic 
short distance potential. The masses and wave functions of the various states are obtained using 
a variational technique, which allows us to compare the results for both perturbative and nonper- 
turbative treatments of the potential. As well as comparing the mass spectra, radiative widths and 
leptonic widths with the available data, we include a discussion of the errors on the parameters 
contained in the potential, the effect of mixing on the leptonic widths, the Lorentz nature of the 
confining potential and the possible cc interpretation of recently discovered charmonium-like states. 



1. INTRODUCTION 



The discovery of the 3/tp in 1974 [l|,0] marked the beginning of an extended interplay between experimental 
results and the phenomcnological treatment of heavy mesons using quantum chromodynamics. The prediction 
and discovery of the various excited states of the charmonium system and later the upsilon system encouraged 
the investigation of a variety of approaches to understand the systematics of quarkonium spectra and decay 
modes. These include potential models [3, H S, @, 0, 0, S, EE E2 El El EI EE El, E3 ^effective field theory 
treatments [3 El, HI SI HI HI and, more recently, lattice gauge theory calculations [H, [H, HI [27[ . While the 
latter promise to provide an accurate nonperturbative description of mesons and baryons, there is still interest 
in perturbative and quasi-perturbative models, particularly in the investigation of the angular momentum and 
spin effects in heavy meson spectroscopy. 

Over the past 25+ years, potential models have proven valuable in analyzing the spectra and characteristics 
of heavy quarkonium systems p8j . Motivation for revisiting the potential model interpretation of the cc and bb 
systems at this time is provided by recent experimental results: 

1. The discovery of several expected states in the charmonium spectrum (r] c (2S) and h c (lP)) 

2. The discovery of new states (X(3872), X(3940), F(3940), y(4260)) which could be a interpreted as above 
threshold charmonium levels 

3. The discovery of the 1 3 D2 state of the upsilon system 

4. The determination of various E\ widths for cc and bb. 

Our objective here is to examine to what extent a semi-relativistic potential model which includes all v 2 /c 2 and 
onedoop QCD corrections can fit the below threshold cc and bb data - both spectra and decay widths - and 
accommodate the new above threshold states. 

One of the earliest and simplest non-relativistic potential models for heavy quarkonia is that of Eichtcn, et 
al. the Cornell model. This model contains a long range linear term to provide quark confinement, a 
short distance Coulombdike term inspired by the zeroth-order QCD interaction and a treatment of states above 
the continuum threshold. In this approach, splittings within triplet state multiplets of a given orbital angular 
momentum L are generated indirectly by differences in the couplings of states with different total angular 
momentum J to the continuum. By including relativistic corrections to the potential corresponding to scalar or 
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vector particle exchange [1, H[ , it is possible to generate interaction terms depending on orbital and spin angular 
momenta of the quarks without introducing additional parameters. These interaction terms provide a direct 
mechanism for splitting the various J states of a triplet angular momentum multiplct. 

Further refinement of the potential requires a consideration of how the orbital angular momentum and spin- 
dependent interactions of the quarks are modified by the qua rk-g luon interaction of QCD. The full one-loop 
quark-antiquark interaction potential was obtained [29l.j30l. l3ll I32T l33l] and employed in calculating the spectra 
and decays for the cc and bb systems in non-relativistic [7j and relativistic model calculations 

HQ. The spectra 

determined in these models proved to be very accurate for both cc and bb. In the case of the T system, the 
addition of the one-loop corrections enabled an accurate prediction of the energies and splittings of the IP, 2P 
and ID levels 0- 

More recently, the precise determination of leptonic and El decay widths, along with the discovery of several 
new states, some expected and some not, has led to a renewed interest in the potential model description of heavy 
quarkonia [H, [H, Q [IE E3, HH Ull • We have revised and extended the approach of our earlier papers in order 
to investigate the newly measured states and decays as well as to discuss other questions of modelling interest. 
These include the scalar/vector mixture of the phcnomcnological confining potential, errors associated with the 
determination of the parameters appearing in the potential, and the accuracy of the perturbation expansion as 
determined by treating the Hamiltonian as "complete" in the variational calculation discussed below. 

In the next Section, we describe the potential model in some detail. This is followed, in Section[3j by an outline 
of our calculational approach. In Section 31 we present our results for the charmonium and upsilon systems, and 
then give some conclusions in Section [5] Some calculational details are given in the appendix. 



2. SEMI- RELATIVISTIC MODEL 



In our analysis, we use a semi-relativistic Hamiltonian of the form 



H = 2Jp 2 +m 2 + Ar~ — 
= H + V L + V S , 



2ir D7T 



Vl + Vs 



(1) 
(2) 



where /i is the renormalization scale, rif is the effective number of light quark flavors and je is Euler's constant. 
Vl contains the v 2 /c 2 corrections to the linear confining potential 
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where A is the linear coupling strength. The first line in Eq.(3) is the contribution from scalar exchange while 
the second line is the contribution from vector exchange, with fy representing the fraction of vector exchange 
in the interaction. The short distance potential is [1^, [3(| Hl[ 



V s = Vhf + V LS + V t + V S i, 



(4) 
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We have chosen such that it contains the relativistic kinetic energy and the leading order spin-independent 
portions of the long-range confining potential and the one-loop QCD short-range potential. Note that, unlike all 
other terms in Eqs. (O, the last term in Eq. (|5dj) has a 1/m behavior. 



CALCULATIONAL APPROACH 



The cc and bb mass spectra and corresponding wave functions are obtained using a variational approach. The 
wave functions were expanded as 



n 
fc=0 



k+e 



(6) 



where 3^ s (^) denotes the orbital-spin wave function for a specific total angular momentum j, orbital angular 
momentum i and total spin s. The C^'s are determined by minimizing 
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(7) 



with respect to variations in these coefficients. This procedure results in a linear eigenvalue equation for the 
Cfc's and the energies (see Eq. (|A1|1) and is equivalent to solving the Schrodinger equation. The wave functions 
corresponding to different eigenvalues are orthogonal and the k th eigenvalue Afc is an upper bound on the exact 
energy Ef.. For n = 14, the lowest four eigenvalues for any £ are stable to one part in 10 6 . 

We performed the calculations in two ways: (i) perturbatively, as implied above with Hq as the unperturbed 
Hamiltonian and all other terms treated as first-order perturbations; (ii) nonperturbativcly, with all terms 
included in the unperturbed Hamiltonian. The most significant effect of the different treatments is on the wave 
functions. Both approaches yield acceptable spectra, but the leptonic decay widths (which depend strongly on 
the behavior of the wave functions at the origin) turn out to be better in the perturbative treatment. On the 
other hand, the non-perturbative treatment provides a better description of certain El decays. The results of 
both approaches are shown below. 

In either approach, we determine an optimal set of potential parameters a = (a\, oli, - • • , a n ) by minimizing 
the x 2 function 
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where the Oi denote the experimental and theoretical values of some quarkonium observable and the Oi are the 
associated errors. In this work, the Ooxpi consist of a subset of the measured onia masses and leptonic widths. 
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For the masses, the <7j are taken to be the actual experimental error and a common intrinsic theoretical error 
added in quadrature. The latter error reflects the theory uncertainty associated with omitting corrections beyond 
one- loop and is estimated by requiring the x 2 /degree of freedom to be approximately unity. Typically, this error 
is a few MeV. The minimization of \ 2 with respect to variations of the parameters a is accomplished using the 
search program STEPIT [H|. 



4. RESULTS 

We summarize our results in the following tables. It should be noted that the parameter sets used in the 
charmonium and upsilon fits are independent of each other. The parameters resulting from our fits are given in 
Table [I] The errors on A, ag and m q are those that change \ 2 by approximately 1 when the parameters other 
than the one in question are allowed to vary. 
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2.60 


1.32 


4.74 
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0.00 
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TABLE I: Fitted Parameters for the cc and bb systems 



As can be seen from Table HI the value of the parameter A characterizing the confinement strength is a smaller 
for the perturbative treatment compared to the non-perturbativc treatment, but the variation is not great. The 
differences between perturbative cc and bb values for A arc rather small, as are the differences between the 
non-perturbative values. This may be an indication of some slight flavor dependence of the confining term. The 
values of as in the two treatments are very nearly identical for a given qq system and the difference between the 
cc and bb values is consistent with the fact that the scale in the upsilon system is higher. In each case, the value 
of the quark mass is smaller for the perturbative treatment compared to the non-perturbative treatment. The 
primary difference between the parameters obtained for the perturbative and non-perturbative approaches is in 
the values of the scale parameter /i and the fraction fv of vector exchange included in the confining potential. 
In both systems, the value of fv is zero in the perturbative case and of order 20% in the non-perturbative case. 
Similarly, the renormalization scale is always larger for the perturbative treatment. We made no effort to restrict 
the value of /x while fitting other than to be sure that it didn't settle on a value that would suggest a change in 
the number of active flavors lighter than the quark flavor being considered. 



4.1. Charmonium 

The results for our determination of the charmonium levels are shown in Table |TT1 where the * denotes the 
states used in the fit [13 ■ We examined the effect of level mixing induced by the tensor interaction, and found 
that it has a very modest effect on the spectrum. However, as discussed below, s-d mixing does have a significant 
effect on leptonic decay widths. 

In both the perturbative and non-perturbative treatments, the overall fit to the spectrum is quite good, 
with the h c (lP) mass well described and the r/ c (2S) mass somewhat low. Neither approach does very well in 
reproducing the ip(3S) level, although the i/j(4:S) level is accurately predicted. 

Recently, a number of states, which could be interpreted as above threshold cc states, have been observed in 
B decays, at e + e~ colliders and at hadron colliders. Among these, the X(3872) is the most firmly established 
[Hjllsi, with m(X(3872)) = 3871.2 ± 0.5 MeV. An analysis of the quantum numbers of this state [40] shows that 
the assignments J PC = 1 ++ , 2 h are the only ones capable of describing the data. If interpreted as charmonium 
states, the assignments are either Xci(2P) or 1 1 Z?2- The former is more consistent with the masses in Table HT1 
given the tendency of the potential model to predict a larger value for above threshold states. 
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TABLE II: Perturbative and nonperturbative results for the cc spectrum are shown. The perturbative fit uses the indicated 
states and the leptonic widths of the i/j(1S) and if)(2S). In the nonperturbative fit the r] c (2S) and i)j(1D) are included 
and no leptonic widths are used. 

The Belle collaboration reported two states, the X(3940), observed [4l[ in the recoil spectrum of e + e~ — > 
J/tpX, and the y(3940), observed [42| in the decay B — > K(ujJ/ip). These appear to be different states |4lj . 
with the former also seen to decay into D* D. Since the other states seen recoiling against the J/i/j have J = 
[4lj . it is tempting to associate J = with the X(3940) [43|. The charmonium assignments could then be r](3S), 
which would imply an unusually large 35" hyperfine splitting, or the x c o(2P), whose predicted mass tends to be 
too low. With the range of x c j(2P) masses in Table HT1 the F(3940) could be accommodated as a 2P state. 

The BABAR collaboration reported a state y(4260) [44[ observed in the initial state radiation (ISR) reaction 
e + e~ — * 7is 7r + 7r _ J ftp and this state was confirmed by the CLEO collaboration [4{|. The ISR reaction assures 
J PC = 1 , which, if the Y(A26Q) is a conventional charmonium state, implies an assignment of ip{nS) or 
ipinD). It has been suggested HU that the F(4260) rather than the ^(4415) be identified with the 0(45"), but 
this is at odds with both the absence of any 1 level in the vicinity of this mass in most potential models and 
the absence of a large 7r + 7r~ J/tp signal from the '0(4040) [45| . 

Included in Tabic [Til is the lowest multiplet of /-states. The center of mass of these states is essentially the 
same for both the perturbative and non-perturbativc treatments, lying somewhat below the i/j(3S). However, 
the pattern of splittings between the 3 Fj levels is quite different in the two cases. 
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In computing the spectrum variationally, the wave functions are sampled in an average sense. The radiative 
and leptonic decay widths are of interest precisely because they are sensitive to more detailed features of the 
wave functions, either matrix elements of the dipole operator in the case of El decays, or the behavior of the 
wave function at the origin in the case of leptonic decays. 

As is usual in potential model treatments 0, [TH, [l6j|, the radiative widths were calculated in the dipole 
approximation. We obtained the E\ and M\ matrix elements by using the variational radial wave functions to 
construct initial and final state wave functions with the appropriate angular dependence and explicitly performing 
the angular integration. Our results are equivalent to the formulas 



T(n 3 Lj - n' 3 L'j, + 7 ) = ^WC /4 | (»' 3 L' r \r\n 3 Lj)\ 2 g , (9) 



for Ei transitions, and 



T(n 2s+1 Lj ^ n ' 2s ' +1 Lj, + 7 ) = V " o, 3 ^-^<W ±1 ^ , s + i Lj/ \ n *>+i Lj) , (i ) 



4 , a ,2J' + 1 



3 to 2 2L + 1 



2 E 



Mi 



for Mi transitions. Here, to is the photon energy, q is the quark charge in units of the proton charge, Ef is the 
energy of the final quarkonium state, Mi is the mass of the initial quarkonium state, m q is the quark mass and 
C fi contains the square of a 6j symbol that arises from the angular integration [47| . 

The resulting cc radiative widths are shown in Table IIIII for both the perturbative and non-perturbative 
treatments. Comparing the two approaches, it is noteworthy that there is a clear indication of a difference in 
behavior of the p state wave functions. This is most easily seen in the results for the ip(2S) — > Xc jO-P) transition 
widths. In the perturbative case, the dipole matrix element is that same for every J, and the widths reflect the 
differences in phase space. On the other hand, the non-perturbative widths are very nearly equal, which means 
that the matrix elements for the various values of J must differ. This difference is also reflected in predictions 
for the ip(lD) — > Xc jO-P) 7 decays recently reported by the CLEO collaboration [i^ . 

Unlike the radiative transitions, which are sensitive to matrix elements of the dipole operator, the leptonic 
widths are sensitive to the behavior of the wave function at the origin. Quite generally, the expression for the 
leptonic width of a l -- state can be written 

m- * 87raV|(0|J(0)|O| 2 

T(l -> ee) = — 5 11 

3 m% 

where (0\J(Q)\ipf I ) is the matrix element of the quark current iq(0) , y fl q(Q) between the vacuum and the 
bound state. If ^ a /3(x) is the bound state two-body wave function, it is shown in Ref. [4!| that 

(0|J"(0)K) = iTr[C-V*(0)], (12) 
where C is the charge conjugation matrix. This leads to the expressions (49| 



T eg (^nS)) = ^H^d 2 TesMnD)) = 1 ^%2V ■ (13) 

m nS ^ m nD 



The integrals I n o and I n 2 are 



InO = \^ r dpp 2 \\ + \ ^ I cpM In2 = r dpp* \ 1 - ] ^(p) . 

VttJo \Z 3 J p 2 +m 2j VttJo ^ Jp 2 + m 2 q 



Here, 4> n o(p) and 4> n 2{p) arc the momentum space wave functions and the factors involving yp 2 + to 2 , represent 
relativistic corrections. To the leading order in p 2 /to 2 ,, these expressions reduce to the familiar forms 



r^MnS)) = -^\ Rn0 (o)\ 2 r eS ^(nD)) = 2 q 4 \K 2 (o)\ 2 ■ (is) 

'nS zm nD m q 



For charmonium, we find that the relativistic effects on the s states are rather small for both the perturbative and 
non-perturbative cases. The d states, on the other hand, receive large corrections in both cases, with the non- 
perturbative case receiving the larger of the two. To be consistent, we have used the expressions in Eq. (|13p for all 
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TABLE III: The the radiative decays of the charmonium system are shown. The tf)(lD) — > xj(1P) widths marked with 
a * are from [iH; see also [43( |. 
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TABLE IV: The leptonic widths of the J = 1 states are shown. 



charmonium leptonic decay widths. For the s states, we also include the QCD correction factor (1 — 16as/37r). 
As mentioned above, we included mixing in our perturbative treatment and found little effect on the spectrum 
due to small mixing angles. However, even a small s-d mixing angle has a noticeable effect on ijj(nD) leptonic 
widths, and we included this mixing in the perturbative calculation. The results are shown in Table HVl The 
message here is rather mixed in the sense the perturbative results for the ip(lS) and ip(2S) widths, which were 
included in the fit, are in reasonable agreement with the data, while the ?p(3S) and ipi^S) widths are too large. 
Despite the inclusion of mixing, the ?/>(l-D) and ip(2D) widths are small compared to the data. On the other 
hand, the non-perturbative treatment, which has no leptonic widths in the fit, does a reasonable job accounting 
for the ip(3S) and V'(45') widths, but predicts i/j(1S) and widths that are too small. With the inclusion of 

the relativistic correction, the non-perturbative ip(lD) and tp(2D) widths are of the right order of magnitude. 



4.2. Upsilon 

The results for our determination of the upsilon levels are shown in Table El where, again, the * denotes the 
states used in the fit. In this case, both approaches give very good fits to the spectrum. In the perturbative 
case, the determination of the parameters is quite robust in the sense that fits using fewer states than the eight 
indicated in Table W\ yield parameters, masses and values of % 2 that are quite similar to the ones listed. For 
instance, using only the T(1S),T(2S),T(3S) and T(45) states results in a fit in which the l 3 Pj multiplct is 
displaced upward by about 7 McV and the rest of the spectrum is well described. The nonperturbative fit 
includes the leptonic width of the T(IS') and the main difference in the results of the two approaches occurs 
in the T(nS) — 7]b(nS) hypcrfinc splitting, which is always smaller in the nonperturbative treatment. In the T 
system there are / states occurring below BB threshold. The l 1 ^, l 3 Fj multiplet, in the vicinity of the T(35), 
is certainly below threshold and the 2 1 F3, 2 3 Fj multiplet is very near threshold. 

Our results for the upsilon radiative transitions are shown in Table IVII Due to a lack of knowledge of the 
widths of the Xbj(l-P) an d Xb j(2P) states, direct comparisons to experimental data are limited primarily to 
T(nS) — > Xbj(, n 'P) transitions. For these, both the perturbative and nonperturbative wave functions are able to 
account for the radiative widths quantitatively, including the 61±23eV width of the T(35) — ► Xbo(lP) transition. 
Regarding the Xbj(nP) -> T(n'S)) decays, one can compare T( X bj(2P) -> T(2S))/T(xbj(2P) -> T(15)) to the 
ratios of the measured branching ratios (T1/T2 in Ref. [l^)- These are appended to the bottom of Table IVT1 
Overall, both approaches give a good description of the data, with the perturbative results getting the nod in 
the case of T(nS) — > Xbj( n 'P) transitions and the nonperturbative results being more consistent with the Ti/r2 
ratios. Given the current errors, though, the results from both approaches are certainly satisfactory. 

The calculated leptonic widths for the T(nS) states, including the QCD correction, are compared with the 
experimental data in Table IVIII For the T system, the relativistic corrections of Eq. (fT3|) are negligible. In the 
perturbative case, none of the leptonic width data were used in the fitting and the agreement is very good. For 
the nonperturbative case, where the T(IS) leptonic width is used in the fit, the results are almost identical. 
Unlike the perturbative case, failure to include the T(1S) leptonic width in the fit leads to an unsatisfactory 
description of all the T(nS) leptonic widths. The calculated leptonic widths of the T(1D) and T(2D) are very 
small 0.02 kcV in both perturbative and nonperturbative treatments. 
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m bb (MeVJ 


r^ert 


Non-pert 






9413.70 


9421.02 




T(15)* 


9460.69 


9460.28 


9460.30 ± 0.26 




9861.12 


9860.43 


9859.44 ± 0.52 




9891.33 


9892.83 


9892.78 ± 0.40 


X62(1P)* 


9911.79 


9910.13 


9912.21 ± 0.40 


/l 6 (lP) 


9899.99 


9899.94 




r?fc(2S) 


9998.69 


10003.6 




T(25')* 


10022.5 


10023.5 


10023.26 ± 0.31 


T(1D) 


10149.5 


10148.8 




1 S D 2 


10157.1 


10157.0 


10161.1 ± 1.7 


I s D a 


10162.9 


10164.1 




I D 2 


10158.4 


10158.3 






10230.5 


10231.4 


10232.5 ± 0.6 


Xbi(2py 


10255.0 


10257.6 


10255.46 ± 0.55 


Xb2{2py 


10271.5 


10271.1 


10268.65 ± 0.55 


h b (2P) 


10262.0 


10263.1 




1 '-i -w — 1 


10353.0 


10351.0 




1 3 F 3 


10355.8 


10355.6 




I 3 Fa 


10357.5 


10359.7 




1ft 


10355.9 


10355.9 






10344.8 


10350.4 




T(3S) 


10363.6 


10365.6 


10355.2 ± 0.5 


T(2D) 


10443.1 


10443.7 




2 3 D 2 


10450.3 


10451.2 




2 3 D 3 


10455.9 


10457.5 




2 1 I> 2 


10451.6 


10452.4 




2 3 F 2 


10610.0 


10609.0 




2 3 F 3 


10613.0 


10613.4 




2 3 _p4 


i nfii ^ n 

1UU1J.U 


i nisi 7 i 




2 : ft 


10613.2 


10613.7 






10622.8 


10631.5 




T(45) 


10643.0 


10643.4 


10579.4 ± 1.2 



TABLE V: Perturbative and nonperturbative results for the 66 spectrum are shown. The perturbative fit uses the indicated 
states . 

5. CONCLUSIONS 

We have shown that a potential model consisting of the relativistic kinetic energy, a linear long-range confining 
potential together with its v 2 /c 2 relativistic corrections, and the full v 2 /c 2 plus one-loop QCD corrected short 
distance potential is capable of providing extremely good fits to the spectra of the cc and bb heavy quarkonium 
systems. Interestingly enough, the the results are about as good with the spin-dependent and non-leading 
order spin-independent terms of the potential treated either as a perturbation or as part of the unperturbed 
Hamiltonian. We find that for both the charmonium and upsilon systems, the perturbative treatment requires 
the long-range potential to be entirely due to scalar exchange, while the non-perturbative treatment requires the 
long-range potential to be about one-fifth vector exchange. 



i -y [Ke V ) 


.rert 


Non-pert 


rjjXpt 


T(1S) -> r?i,(15) 


0.004 


0.001 




T(25) -» ?7i>(2S') 


0.0005 


0.0002 




T(25) -» ^(lS 1 ) 


0.0 


0.005 


< 0.02 


T(25) -» X6o(lP) 


1.15 


0.74 


1.22 ± 0.16 


T(25) -»X6i(lP) 


1.87 


1.40 


2.21 ± 0.22 


T(25) ^X62(1P) 


1.88 


1.67 


2.29 ± 0.22 


%(2S) -> iifc(lP) 


4.17 


20.4 




T(3S) -> Xfto(2P) 


1.67 


1.07 


1.20 ± 0.16 


T(35) Xi>i(2P) 


2.74 


2.05 


2.56 ± 0.34 


T(35) -> X6 2(2P) 


2.80 


2.51 


2.66 ± 0.41 


T(35) -> Xbo(l-P) 


0.03 


0.03 


0.061 ±0.023 


T(35) -> xm(I-P) 


0.09 


0.003 




T(3S) -> Xf>a(lP) 


0.13 


0.11 




Xi,o(lf) ->• T(15) 


22.1 


19.6 




Xfci(l^) T(15) 


27.3 


23.9 




Xb2(l-P) T(15) 


31.2 


26.3 




ftb(lP) -> 7?6 (15) 


37.9 


4.61 




X6o(2P) -> T(25) 


9.90 


9.91 




X6i(2P) -> T(25) 


13.7 


12.4 




X6 2(2P) -> T(25) 


16.8 


13.5 




X6o(2P) T(15) 


6.69 


1.83 




X6i(2P) T(1S) 


7.31 


4.81 




X6 2 (2P) -> T(IS') 


7.74 


6.86 




X6o(2P) -> T(1.D) 


1.13 


1.05 




X6i(2P) T(LD) 


0.62 


0.52 




X6 2 (2P) T(LD) 


0.04 


0.03 




X6i(2P) -> 1 3 L> 2 


1.48 


1.31 




X62(2P) -> 1 3 £>2 


0.47 


0.35 




T(1D) ^Xbo(lP) 


18.1 


12.5 




T(LD) ^X6i(lP) 


9.82 


7.59 




T (1Z?) X62(1P) 


0.51 


0.44 




± 3 D 2 -> Y6l(lP) 


19.3 


14.9 




1 3 D2 -» YfcaflP") 


5.07 


4.35 




1 3 7J 3 _> y h2 (lP) 


21.7 


18.8 




Ti/r 2 


Pert 


Non-pert 


Expt 


ri(x6o)/r 2 (xbo) 


1.48 


5.42 


5.11 ±4.14 


ri(x6i)/r 2 (x6i) 


1.87 


2.58 


2.47 ±0.60 


ri(xb2)/r 2 (xb2) 


2.17 


1.97 


2.28 ± 0.47 



TABLE VI: The radiative decays of the upsilon system are shown. 
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1 e e [KeV ) 


.r ere 


Non-pert 


rjxpt 


T(1S) 


1.33 


1.33 


1.340 ± O.Olo 


T(2S) 


0.61 


0.62 


0.612 ±0.011 


T(3S) 


0.46 


0.45 


0.443 ± 0.008 


T(4S) 


0.35 


0.30 


0.272 ± 0.029 


TABLE VII: 


The leptonic 


widths of the T(nS) states 


are shown. 



The photon and leptonic widths obtained from the variational wave functions are, for the most part, in very 
good agreement with the available data. An interesting difference between the two treatments is that in the 
non-perturbative treatment every state has its own wave function, while in the perturbative treatment all states 
in the same angular momentum multiplet have the same wave function. It is this feature that is responsible 
nonperturbative treatment's somewhat better description of the radiative widths. In the charmonium system, 
the leptonic widths are better described by the perturbative wave functions, while in the upsilon system, the 
nonperturbative description is just as good as the perturbative description provided one leptonic width is included 
in the fit. In both the charmonium and upsilon systems, inclusion of the QCD correction to the s state leptonic 
widths is essential to obtain agreement with experiment. 

Although in some respects the perturbative treatment, and hence pure scalar exchange for the long range 
potential, yields a somewhat better description of the charmonium system, the non-perturbative approach, 
with mixed scalar and vector exchange, is also viable. In the upsilon system, both approaches give very good 
descriptions of the available data and, given the rich spectrum below the continuum threshold, it should be 
possible to decide which is preferable with additional data such as the T(nS) — rjb(nS) hyperfine splittings. 
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APPENDIX A: SOME CALCULATIONAL DETAILS 

The choice of the variational wave function, Eq. $Q results in an equation for the CVs and the energies of the 
form 

n n 

Y,(i\H\k)C k = \J2(i\N\k)C k , (Al) 

k=0 j=0 

where for a state with fixed j,£ and s 

(i\H\k) = Jd 3 r(^y + \- r / R (y^(n)yH^y + \- r / R y^(n), (A2) 

// r \ i+k+2£ , + K>3 

d 3 r(-) e -2r/« (y™ (O)) f Jff,(n) = 2 , +fc+2 , +3 r(i + k + 21 + 3) . (A3) 

The matrix elements (i \ V\k) of the potential, Eqs. ([3]) and (|4]), can be evaluated analytically in coordinate space. 
The contribution of the kinetic energy, 2y / p 2 + m 2 , is evaluated numerically in momentum space, where the 
expansion of the wave function is 

= (~i) l R 3 E c * J n9 ^or(k + 21 + 3)p- + f 7 + 1 3 / / 2 2 (cos 0)y;g B (n) . (A4) 

fc=0 
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Here, P k + e +0 2 ( cos ®) ^ s the associated Legendrc function and 9 is related to the magnitude of p as 

tan 8 = pR . 



(A5) 



Given the matrix elements (i \H\k) and (i \N\k), eigenvalues and eigenvectors for the linear system Eq. (|A1[) 
can be obtained using any standard package [50] ■ We typically use 14 terms in the expansion Eq. ([6]). 

Calculation of the matrix elements of potential using the variational wave function Eq. ^ is relatively straight 
forward. Perhaps the only integrals requiring special attention are those involving 



/ drr 2 \7 2 


In fir + 7_e 


10 


r 



F(r), 



(A6) 



where F(r) denotes a product of the radial wave functions, which has the form (r y 'R) L+L ' e~ 2r '< IR . UL + L'^0, 
then, effectively, 



In /j,r + je 



(A7) 



and the integral in Eq. (|A6j) reduces to a Gamma function. When L + L' = 0, we have 



1 

4* Jo 



drr 2 \7 2 



In fir + ~f E 



- 2r ' R = 1 - In 



~~2~ 



(A8) 



When calculating the matrix elements of the 8(f) terms in Eqs. (|5aj) and (|5d|) , we 'soften' their singularity by 
adopting the quasistatic approximation of Ref. [lj| , which leads to the replacement 



5(f) 



(A9) 



where m is the quark mass. This softening helps the stability of the eigenvalue calculation, particularly in the 
non-perturbative approach. 
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